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ABSTRACT 

We use a simple parameterization of the mid-IR spectra of a wide range of 
O \ galaxy types in order to predict their distribution in the Infrared Array Cam- 

^ \ era (IRAC) 3.6 /zm, 4.5 /im, 5.8 /zm, and 8.0 /im and Multiband Photometer for 

■ Spitzer (MIPS) 24 /xm color-color diagrams. We distinguish three basic spectral 

types by the energetically dominant component in the 3-12 /im regime: stellar- 
^ . dominated; PAH-dominated; and continuum-dominated. We have used a Markov 

chain Monte Carlo (MCMC) approach to arrive at a more systematic and robust 
representation of the mid-IR spectra of galaxies than more traditional approaches. 
We find that IRAC color-color plots are well suited to distinguishing the above 
spectral types, while the addition of 24 /im data allows us to suggest practical 3- 
color cuts which preferentially select higher redshift sources of specific type. We 
compare our simulations with the color-color plot obtained by the Spitzer First 
Look Survey (FLS) and find reasonable agreement. Lastly, we discuss other 
applications as well as future directions for this work. 

Subject headings: galaxies: infrared, modeling 
1. Introduction 



^Visiting graduate student at the Spitzer Sci- 
ence Center for the period during which most of The 3-12 /im mid-IR spectral band 
this work was done; sajina@astro.ubc.ca lies at the interface between the regimes 
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where one observes stellar light directly 
and where one sees that light reprocessed 
by dust and gas. The main contribu- 
tors to galactic emission at mid-IR wave- 
lengths are light from old stellar popu- 
lations (Bosselli et al. 2003), Polycychc 
Aromatic Hydrocarbons (PAH) emission 
features (Puget et al. 1985, Hudgins & 
AUamandola 2004), and power law con- 
tinuum emission, which is, in particular, 
the dominant component in the mid-IR 
spectra of type I AGN mid-IR spectra 
(Clavel et al. 2000). This last constribution 
is likely due to stochastically-heated very 
small grains (Desert et al. 1990). These 
sources of emission can be extincted by 
dust, and further attenuated by discrete 
absorption features, in particular those 
from water ice (at ~6/xm), and Si (at 
~10/im). Other possible contributors in- 
clude: Asymptotic Giant Branch (AGE) 
stars' dust shells; a near sublimation tem- 
perature (~ 1000 K) blackbody component 
(Lu et al. 2003); molecular cooling lines; 
and nebular hues. Given this complex 
spectrum, how much Spitzer broadband 
photometry alone can tell us about the 
nature of the sources, especially if the red- 
shift is also unknown? And how can one 
interpret the IRAC (Fazio et al. 2004) 
color-color diagram (specifically S^^g/Ss.e 
vs. 5'8.o/'S'4.5), as well as related diagrams 
including the MIPS (Rieke et al. 2004) 
24 fim band? 

One approach to the above problem 
would be to use a model with a large num- 
ber of parameters accounting for all known 
sources of emission and attenuation; how- 
ever, this complexity is not justifiable given 
the small number of available photometric 
data points. Another approach would be 



to use a library of fixed spectral templates 
to choose between, which will allow for the 
general trends to be studied. However, lit- 
tle information can be extracted from the 
distribution of observed colors unless a well 
defined set of templates is used. 

We choose a different approach where 
we first search for a common empirical 
description of the ~ 1-12 fim spectral en- 
ergy distributions (SEDs) of galaxies, fo- 
cusing only on what drives the basic shape 
rather than the detailed structure. In 
our model, the SED of essentially any 
galaxy can be representated by a combina- 
tion of direct stellar emission, PAH emis- 
sion, and a power law continuum. We fit 
archival ISO spectra of a range of galaxy 
types with this model using Markov Chain 
Monte Carlo (MCMC) approach, which 
provides a handle on the posterior prob- 
ability distribution of each parameter. In 
addition, the Markov chains resulting from 
fitting the spectra of a wide range of galaxy 
types represent a distribution of poten- 
tial spectra consistent with the input ones. 
We simulate the IRAC and IRAC+MIPS 
color-color diagrams including redshift ef- 
fects (2; ~ — 2) by drawing from the gen- 
erated probability distributions, while in- 
cluding an additional noise term to com- 
pensate for the simplicity of the model. 

Lastly, we compare our predictions with 
early results from the Spitzer Space Tele- 
scope, and discuss future directions for this 
work. 

2. The sample 

In order to study how a wide range of 
galaxy types map onto Spitzer photome- 
try bands, we need to use a local sample 
for which the mid-IR SEDs are determined 
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with sufficiently high spectral resolution 
and signal-to-noise ratio (SNR). For this 
purpose, we use a sample of 60 local galax- 
ies, including ellipticals, normal spirals, 
ultraluminous infrared galaxies (ULIGs), 
starbursts, Seyfert galaxies, and HII galax- 
ies (see Table 1). This sample is in no 
sense complete or free of selection bias (the 
most obvious of which is toward more ac- 
tive galaxies), our approach being merely 
to cover as wide a range of spectra as pos- 
sible. The limiting factor was the avail- 
ability of 2-12 fim spectroscopy, especially 
for classes of galaxy whose volume den- 
sity in the local Universe is low. We use 
spectra taken with the PHT-S spectropho- 
tometer (2.5-4.9 fim and 5.8-11.6 /xm, with 
A/AA=100-300) on ISO, which has a 24" x 
24" aperture. For comparison, IRAS had 
resolution from about 30" at 12 //m to 
about 120" at 100 //m. We refer to Lu 
et al. (2003) for details of the ISO spec- 
tra and data reduction. We only include 
in our sample sources detected in all 4 
IRAS bands, which naturally restricts us to 
bright, and/or local galaxies ((2;)=0.009). 
The typical SNR is ~200\ however, some 
of our sources have significantly higher 
SNRs which is likely to bias the fitting pro- 
cedure and certainly makes it difficult to 
interpret the best achieved. To avoid 
this bias, we have reweightcd the spectra 
with weights corresponding to a signal-to- 
noise of 200. 

The above spectra were extended into 
the near-IR (NIR) by extracting the flux 
from 2MASS images using an aperture 
equivalent to the PHT-S one in a manner 
similar to that of Lu et al. (2003) ^ 

lAs estimated by E(a;i/^|)/VE(1M')- 

^The main difference being that we adopted the 



In 10-20% of the sample the last ~ 0.4 /xm 
of each half of the PHT-S spectra (i.e. 4.5- 
4.9 /im and 11.2-11.6 fim) suffer from what 
appear to be calibration problems, thus for 
simplicity we ignore these sections in the 
subsequent analysis. All spectra have been 
/c-corrected to redshift zero. 

3. Model 

As stated in the introduction, the ma- 
jor components governing the shape of the 
1-12 /xm SEDs of galaxies are: direct light 
from old stellar populations; a smooth con- 
tinuum probably due to very small grains; 
and PAH emission features. Fig. 1 illus- 
trates this idea. Our model treats all other 
possible contributors as noise to the above 
three. We now discuss each of the three 
model components in turn. 

Stars: We assume a 10 Gyr old, Salpeter 
IMF, single stellar population (SSP) which 
is justifiable since we only consider the 
SED redward of 1 /^m where it is domi- 
nated by K and M giants. At the red- 
shift range we consider, we are likely dom- 
inated by evolved systems, and therefore 
assuming solar mctallicity is reasonable. 
Lower- met alhcities would lead to a less 
pronounced 1.6 /xm bump and 2/xm break 
(for a discussion of how these appear to 
IRAC function of redshift see Simp- 
son & Eisenhardt 1999). We consider 
the effect of including the dust shell emis- 
sion of AGB stars in Section 6.1. We 
use the PEGASE2.0 code (Fioc & Roca- 
Volmerange 1997) to generate a spectral 

2MASS flux calibration from Cohen et al. (2003) 
without additional offsets. Where available, the 
images used were taken from the 2MASS Large 
Galaxy Atlas (Jarrett et al. 2003), and the rest 
from the 2MASS Extended Source Catalog. 
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template with the above specifications. 

Continuum: The continuum is given by 
a flux density fiy oc z/~". Note that, un- 
hke here, a power-law is often used to de- 
scribe the entire MIR spectrum. It is likely 
that for the majority of galaxies the ori- 
gin of this component is the non-thermal 
emission of very small grains (Desert et 
al. 1990), although the Wicn tail of a 
> 100 K thermal component is also a pos- 
sible contributor. Lastly, in some extreme 
AGN-dominated sources, it can be due to 
synchrotron emission, although dust tori 
are a more likely source of MIR emission 
in typical AGN (Pier & Krolik 1993). 

PAH features: The PAH bands are 
modeled as a set of Lorentz profiles^ 

(Boulanger, Boissel, Cesarsky, & Ry- 
ter 1998) at 3.3, 6.3, 7.7, 8.6 and 11.3 fim. 
We assume that typically the gap between 
the 8.6 and 11.3 /im PAH features accounts 
fully for the 10 //m gap, without the need 
to invoke Si absorption as well (Hclou et 
al. 2000). For the most active galaxies, the 
usual 7.7/im and 8.6 ;um features can be 
replaced by a broad peak at ~ 8 /xm. This 
may be merely continuum bracketed by 
the water ice and Si absorption features 
(Spoon et al. 2002). However, the posi- 
tions, shapes, and relative strengths of the 
PAH features can vary in different (partic- 
ularly extreme) environments and obscu- 
ration levels. Thus the PAH-dominated 
nature of this feature cannot be ruled out 
(Peeters et al. 2002). For simplicity, we 
assume all such features are due to PAH 
emission, although it is likely that both 
effects play a role. Detailed modeling of 
the PAH features is beyond the scope of 

Ke. f.oc{na{l + {{u-uo)/<Tf)}-' 



this work, the template used being merely 
a fiducial which does not lead to obvious 
'features' in the residuals of the fits (i.e. it 
is some average of the PAH emission fea- 
tures for our sample). 

The overall model can be expressed as: 

= ai X Fi (stars) + as x F2(PAH) + 
10^^ X i>~^^, where the a's designate the free 
parameters, and the F's are spectral tem- 
plates. The last two parameters are used 
to describe the power law component. 

4. Method 

We use a Markov Chain Monte Carlo 
(MCMC) approach both to find the best- 
fit set of parameters for our data-set and 
in order to sample the posterior parameter 
distributions. We describe only the most 
relevant details here, but see e.g. Panter 
et al. (2003) and Lewis & Bridle (2002) for 
more details of the procedure as well as ex- 
amples of other applications in astronomy. 
Starting with uniform priors in our 4 pa- 
rameters, subsequent trial sets are drawn 
from a Gaussian around the last accepted 
set (the proposal distribution, q'^). We fol- 
low the Metropolis-Hasting acceptance cri- 
terion: tti+i = mm[u,pi+iqi+i/piqi], where 
the stochastic element is provided by the 
random number Me[0,l]. Since the pro- 
posal distribution is symmetric, q-i+i/qi is 
unity. The probability (p) that the system 
finds itself in a given state is given by the 
Boltzman factor where the 'energy' is 
here, and thus pj+i/pi is exp(— Ax^/T). 
For a proposal distribution width of 1% of 
the total parameter range, T = 100 results 

■^Note that the exact shape of the proposal distri- 
bution is largely irrelevant, although its 'width' in 
particular affects the efficiency of sampling. 
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in > 50% initial acceptance rate, which 
constitutes good mixing. We use simu- 
lated annealing, i.e. gradually decrease 
T, to allow us to find the high-density re- 
gions more efficiently. We run the Monte 
Carlo for 2000 iterations at each tem- 
perature, which meanwhile is decreased 
as Ti+i = 0.9Ti (Kirkpatrick 1984) for 20 
steps^. Note that the use of simulated an- 
nealing in conjunction with MCMC is non- 
standard. It has the advantage of vastly 
speeding up the procedure, and should not 
severely affect the distribution achieved 
provided the annealing is slow. Although 
we have not tested the reliability exhaus- 
tively, any deviations will be in the tails of 
the distribution and therefore irrelevant for 
the present purposes. There is no univer- 
sal recipe for deciding when convergence 
(the 'burn-in') is achieved (e.g. Panter et 
al. 2003). However, by using a few test 
cases we find that although it is usually 
achieved fairly quickly it depends strongly 
on the initial conditions. We conserva- 
tively cut out the first three temperature 
steps (6000 iterations), which guarantees 
us convergence in all cases. We then thin 
the chain by a factor of 30 to decrease 
the internal correlations which result from 
every link knowing about its predecessor 
(Lewis & Bridle 2002). The result of the 
MCMC analysis is a chain of parameters 
(for each of our 60 galaxy SEDs) which 
properly samples the distribution of best- 
fit model parameters. 

^Thcsc choices are not unique or optimal. We have 
merely chosen a set of MCMC parameters which 
work well for this particular case. 



5. Results 

5.1. Individual fits 

Fig. 2 shows the spectra and fits for a 
sub-sample of our galaxies along with the 
cross-correlation between the stellar and 
PAH amplitude parameters (ai and 02). 
Clearly the model provides a reasonable fit 
to the general shape of the SEDs, although, 
as expected, residual spectral substructure 
remains, largely due to variations in the 
PAH features. This is also seen in the mean 
of the reduced x^s obtained which is ~ 4. 

For our purposes, we need to estimate 
how much noise to include in the Spitzer 
broadband fiuxes in order to account for 
spectral substructure not included in the 
model without including the uncertainties 
associated with the PHT-S spectra them- 
selves. In Fig. 3, we show the histograms of 
the difference in IRAC fiuxes obtained by 
convolving the residuals with the IRAC fil- 
ter profiles as AF — J ASi,Ti,dv/ f T^^du, 
where Tj, is the transmission curve for the 
given filter, and AS^, is the residual spec- 
trum. These distributions peak near zero 
showing that our model accounts fairly 
well for the relevant spectral gradients. 
We estimate the scatter expected from the 
uncertainties in the data by taking their 
average under the filter profiles (analo- 
gously to the flux estimates) and dividing 
by y/N (i.e. the error-in-mean). We com- 
pare this scatter with the one obtained and 
find that the excess is accounted for by 
CTa.e = 10 mJy, 0-4.5 = 8 mJy, (75.8 = 18 mJy, 
and CTg.o = 9 mJy. The higher value for (T5.8 
is due to the lower average uncertainty in 
the data, which is probably related to this 
filter being at the edge of one of the PHT-S 
spectral channels. The observed scatter is 
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similar to (rather than twice) that of the 
other three bands. For this reason, as well 
as because using the value obtained results 
in an unreasonable stretch in the relevant 
colors, we use (Js.g = 10 mJy in the color- 
color plot simulations (Section 6), resulting 
in uniform scatter across the colors. 

5.2. Correlations 

The Markov chains constructed above 
allow us to study how the parameters re- 
late to each other. This is shown in Fig. 2 
where, as expected, the general trend is for 
strong stellar light emitters to be relatively 
weak PAH emitters. 

Fig. 4 shows the marginalized Markov 
chain distributions of the continuum slope, 
a for a representative sub-sample of six 
galaxies. In addition, the total distri- 
butions of the PAH-dominated, stellar- 
dominated, and continuum-dominated sources 
are shown, where the highest peak is, re- 
spectively, a ~ 1.8, a ~ 1.6, and a ~ 1.2 
for the three types. The last is roughly 
consistent with the a = 0.84 ± 0.24 ob- 
tained by Clavel et al. (2000) for Seyfert 1 
galaxies. The continuum component suf- 
fers from some degeneracy with the PAH 
emission strength, making it difficult to 
assess what fraction of the width of these 
distributions is real. 

For the best-fit models for each of our 
sample galaxies, we look at how this pa- 
rameterization relates to the IRAS col- 
ors. Fig. 5a shows that the strength 
of the PAH emission is positively cor- 
related with the Seo/Sioo color, as ex- 
pected since the later is a measure of the 
interstellar radiation field (ISRF) which 
excites the PAH molecules. The excep- 
tions here are the star symbols, which are 



star-formation-dominatcd active galaxies 
(Arp220, NGC4418, 1C860), and the cir- 
cles, which are continuum-dominated AGN 
(NGC4507, NGC3783, Mrk231). This dis- 
tinction is seen in Fig. 5c where the circles 
are clearly strong continuum emitters (al- 
though not the only ones) while the stars 
are weaker than average continuum emit- 
ters (NGC4418 has no discernible contin- 
uum). In Fig. 5d we show the well known 
IRAS color-color plot, and notice that as 
expected the 'stars' occupy the area typ- 
ically explained by very active galaxies 
whose PAH carriers are being destroyed 
by the extreme ISRF. On the other hand, 
the AGN sources obey the general trend, 
showing that these colors are insufficient in 
themselves to distinguish AGN-dominated 
galaxies, as the lack of PAH emission is 
compensated by warmer dust resulting in 
similar S12/S25 colors (IRAC colors are far 
more efficient here, as discussed in Sec- 
tion 7.1). 

Excepting the six galaxies discussed 
above, the ratio of the PAH to stellar emis- 
sion appears anti-correlated with the IRAS 
S12/S25 color (Fig. 5c). This cannot be 
explained through PAH emission, since it 
affects the 12 fim filter (which includes the 
11.3 /xm and 12.6/xm features), and thus, 
if this was the primary driver, the correla- 
tion would be positive. It rather suggests 
that warm (~ 100 K) dust which primar- 
ily affects the 25 fim channel, goes along 
with an increase in the PAH-to-stellar ra- 
tio. This makes sense, since both PAH 
emission and warm dust are expected to 
trace relatively recent star formation (see 
e.g. Forster Schreiber et al. 2004) due to 
their preferentially being excited by ener- 
getic UV/optical photons, while the stellar 



6 



component here is composed of older, red- 
der stars. 

6. IRAC color-color simulation 

Since we have 60 independent spec- 
tra, then the overall probability of the 
data given the model is P{9\D,M) ~ 
nfi^P{9\'Di,M), and we simply link the 
individual chains obtained above to form 
a master chain. Our simulation uniformly 
draws 10,000 times from this chain. It 
is thus subject to the same selection ef- 
fects as our original sample. Each set of 
parameters thus obtained defines an SED 
through the four-parameter model, at a 
resolution of X/6X ~ 200. To obtain the 
IRAC flux densities from these, we con- 
volve the above with the IRAC filter pro- 
files. To these flux densities we add Gaus- 
sian noise, of width equal to the rms for 
the residuals found in section 5.1, in order 
to account for the missing substructure. 
Note that this scatter should be regarded 
as an upper limit as the Markov chain al- 
ready accounts for some of it and therefore 
some fraction of the noise is added twice. 
The resulting IRAC color-color diagram 
is shown in Fig. 6. Here we distinguish 
stellar-dominated (red), PAH-dominated 
(green), and continuum-dominated (blue) 
sources. 

The most prominent feature of this 
color-color plot is the PAH-to-stars se- 
quence. Since, as discussed previously, this 
ratio is an indicator of the star-formation 
activity of a galaxy, this is likely to prove 
a valuable diagnostic in the study of lo- 
cal sources such as resolved galaxies in 
the nearby Universe or even galactic star- 
forming regions. Its exact shape is sen- 
sitive to the PAH feature profile and thus 



can itself be used as a diagnostic for distin- 
guishing for example whether this sequence 
is significantly different in low metallicity 
environments than high metallicity ones. 
As the photometric uncertainty of Spitzer 
observations is typically far smaller than 
the scatter considered here, such studies 
are feasible. 

The second aspect is that continuum- 
dominated sources (AGN) form another 
sequence with redder S^^s/^s.e colors, for 
a corresponding Ss.o/Si,^ color, than the 
above. 

Combinations of our three components 
occupy the region in between these two se- 
quences. Samples with different selection 
criteria than ours will form the same gen- 
eral shape here, but the density of points 
across it will differ. 

The use of mid-IR color-color plots to 
distinguish star-formation-powered from 
AGN-powered sources was already known 
for the ISO colors (Laurent et al. 2000). 
The main difference here is the greater sen- 
sitivity to direct stellar light in the IRAC 
color-color plot. 

6.1. Effect of AGB dust shells 

Assessing the level of contribution of 
AGB dust shell emission to galaxy spec- 
tra is difficult because the effect has the 
same sense (in terms of the IRAC col- 
ors) as interstellar dust. To show this we 
overplot solar metallicity^ lOOMyr-lOGyr 
SSP spectra of Piovan et al. (2003) onto 
our z — color-color plot (Fig. 6). AGBs 
have a significant effect on the spectra 

^Metallicity has a much smaller effect than age 

here, since wc arc mostly sensitive to the relative 
fraction of the SSP of stars in their AGB phase. 
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only for relatively young stellar popula- 
tions. Such populations presumably also 
emit strongly in their ISM dust ( 
sumed in our model) which starts to dom- 
inate the spectra at A > 4 //m. Moreover, 
in spatially integrated galaxy spectra, the 
contribution of the young (e.g. starburst) 
population can be overpowered in the ob- 
served integrated spectra by an older sur- 
rounding population. In Fig. 7 we show 
the actual stellar population SEDs, includ- 
ing AGB dust shell emission and some of 
the available ISO spectra for comparison. 
We notice that, as expected, the difference 
is most easily assessed in the roughly 3- 
5 fjim regime. AGB-dominated SEDs over- 
predict the power emitted in those wave- 
lengths compared with the near-lR emis- 
sion. This suggests that the mean age of 
the stars in the sampled aperture is greater 
than ~ 0.5 Gyr. Some of our power-law 
component may be contributed to by AGB 
emission in the > 5 /zm spectra (see for ex- 
ample 300Myr-old SSP spectra in Fig. 7). 

For the study of integrated galaxy spec- 
tra this effect is not significant. However, 
it should be kept in mind when looking at 
finer physical scales inside galaxies, where 
younger average ages may be observed. In 
addition, possible foreground from Galac- 
tic AGB stars may be an issue specifically 
near the age track in Fig. 6. 

6.2. The effect of redshift 

We show in Fig. 8 the z = 0-2 redshift 
evolution in the IRAC color-color plot for 
models representative of the three basic 
types discussed in this work. We can un- 
derstand the general trends of movement in 
the color-color diagrams through consider- 
ing the main PAH features moving through 



the IRAC bands. 

For PAH-dominated sources, the S's.s/S's.e 
color shifts blueward between 0-0.1, 
due to the 6.2 /xm feature leaving the 
5.8 /xm filter. By 2; ~ 0.2 the 3.3 /xm fea- 
ture leaves the 3.6 /im filter and thus the 
color shifts slightly redward again, until 
the stellar 1.6/im bump begins to enter 
the 3.6;um filter past z ~ 0.7. The color 
shifts redward strongly when the 5.8 /xni 
filters begins to sample the stellar peak 
as well, which makes this color a sen- 
sitive redshift estimator past 2; ~ 1 (Le 
Floc'h et al. 2004). On the other hand, 
the 5'8.o/>S'4.5 color starts ofi' red, due to 
the powerful 7.7//m PAH feature complex, 
and gradually turns bluer, due to the in- 
fluence of the 3.3 fim feature on the 4.5 fim 
filter {z ~ 0.2 — 0.5) coupled with the 
6.2 /xm feature replacing the 7.7/im one. 
Past z ~ 0.5 the strongest PAH features 
leave the 8 //m filter altogether, and at this 
point star-forming galaxies enter the gen- 
eral vicinity of 2; ~ stellar-dominated 
sources, which is also where contamina- 
tion from foreground stars (especially dust 
enshrouded ones) is most severe. Thus, 
these two colors alone are insufficient to 
distinguish the high-2; starforming galax- 
ies. 

Stellar-dominated sources occupy the 
same general vicinity until z ~ 0.6 when 
the 1.6 yum bump enters the 3.6 /im filter, 
while further redshifting (especially past 
z ~ 2) causes a red S5,s/S3s color due to 
their red optical continua probably lead- 
ing to a sequence similar to the 'AGN- 
sequence', however offset in its Ss.o/S^,^ 
color. 

Continuum-dominated sources, being 
relatively featureless, also do not leave the 
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general vicinity they occupy locally, in- 
stead they just shift along the 'continuum- 
sequence' due to slight changes in slope. 
Only largely unobscured AGN will ever 
display 'blue' IRAC colors when the rest- 
frame optical wavelengths begin to be sam- 
pled at 2; > 2. 

6.3. The use of 24 jdui data at 2; > 1 

Due to the lack of available data, we do 
not formally extend our model past 12 /xm. 
However, as seen in the previous section, 
the diagnostic power of the 8 /xm PAH fea- 
ture complex becomes inaccessible to the 
IRAC colors alone past 2;~0.7. We over- 
come this by the use of the MIPS 24 /xm fil- 
ter once our 1-12 /xm spectra are redshifted 
into it at 2;^ 1. As discussed previously, 
the narrow PAH-to-stars sequence is likely 
to be an important diagnostic arising from 
the IRAC color-color plot. Here we con- 
sider the lg(5'24/S'5.8) color which at z > 1 
mimicks the lg(S'8/S'4.5) local PAH-to-stars 
sequence. In Fig. 9, we show the redshift 
dependence of the three colors for the dif- 
ferent spectral types. Obviously, the use of 
such an implicit '24 /xm-detection' criterion 
is also useful in removing foreground con- 
tamination, which is an issue in particular 
as the IRAC colors of high- 2; galaxies can 
be very 'star'-like (see Section 6.2). 

6.4. Color selection 

In Fig. 10, we present specific color- 
cuts in the three colors discussed (S's.g/^'s.e, 
>S'8.o/'S'/i.5, and 5'24.o/5'5.8), which preferen- 
tially select sources with specific types of 
mid-IR spectra at various redshift ranges. 
As in the local IRAC color-color plot, the 
most prominent feature here is the PAH- 
to-stars sequence, which shifts left with 



redshift, allowing both for a handle on 
the underlying A^(2;) and some sense of 
the physical characteristics of the sample 
(since the above sequence is most likely a 
tracer of the SF activity of a galaxy) . Since 
the various regimes discussed can overlap 
significantly, this is likely to be most use- 
ful for large samples where 'mis-hits' are 
washed out, while for small samples (order 
of 10) such color-cuts are hkely to be hard 
to interpret unless additional information 
is available. 

7. Discussion 

7.1. Compcirison with Spitzer re- 
sults 

The First Look Survey (FLS), con- 
ducted with the IRAC and MIPS instru- 
ments on the Spitzer Space Telescope, is 
a pubhcly available, shallow, 4deg^ sur- 
vey^. The IRAC 4-band catalog (Lacy et 
al. 2004b) has fiux density limits (5 a in 
a 5 " aperture) of 0-3.6 ~ 7 fiJj, (74.5 ~ 8 /xJy, 
(T5.8 ~ 60 iiJy, and cs.o ~ 50 /xJy. We use 
this catalog to explore qualitatively how 
our model compares with real data. Fig. 11 
shows the FLS data on the top, with a 
simulation of it from our model shown 
below. To create the latter we use the 
master chain discussed previously. The 
A'"(2;) used covers 2; = 0.2 — 2.0 and is 
a combination of two Gaussians centered 
at 2; = 0.2 and 2; = 0.6 of widths 0.08 and 
0.3 respectively. Note that this is merely 
a curve giving a reasonable approxima- 
tion to the observed FLS IRAC color-color 
plot. Work is in progress to more robustly 
derive the redshift distributions per lumi- 
nosity per spectral type by comparing our 

^http: / / ssc.spitzer.caltech.edu/fls/ 
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simulations with observed color-color dis- 
tributions (see following section). In ad- 
dition, since stars have not been removed 
from the FLS catalog, we add a Gaussian 
distribution centered at (-0.4, -0.6) with 
a— 0.15 and the same centered at (-0.3, 
-0.5) to cover the gap from the above to 
the old stellar population colors used in 
our model, and thus emulate the expected 
continuum of colors from different types of 
stars (no attempt has been made to opti- 
mize this). Although not perfect in detail, 
this nonetheless shows that the essential 
elements of the observed FLS color-color 
plot are well reproduced by our model. 

Our model suggests, that the two 'bunny 
ears' in the FLS color-color plot separate 
the z ~ 0—1 continuum-dominated sources 
(right ear) from the z ~ 0.1 — 0.3 PAH- 
dominated sources (left ear). This ability 
to separate out the continuum-sources (i.e. 
AGN) based on Spitzer photometry, which 
is far less affected by extinction than op- 
tical observations, was explored by Lacy 
et al. (2004a) in order to assess the frac- 
tion of AGN-dominated sources missed by 
optical surveys due to dust obscuration. 

We also want to compare our simula- 
tions with the properties of known z > 1 
sources observed with IRAC and MIPS. 
Ivison et al. (2004) present the Spitzer col- 
ors of nine MAMBO sources. The 1200 //m 
MAMBO filter, by virtue of its negative 
/c-correction (climbing the Rayleigh-Jeans 
tail of the thermal dust emission peak), 
selects higher redshift sources. Ivison et 
al. (2004) indicate how the Spitzer colors 
can be used to separate out AGN from 
star-formation(PAH)-dominated sources. 
We found that, although our sample on 
average did not require extinction applied 



to it to achieve acceptable fits^, this was 
essential to explain the above observa- 
tions. In Fig. 12 we show the Ivison et 
al. (2004) color plot with extinction grad- 
ually applied (using the Draine et al. 2003 
Milky Way Ry — ^-l extinction curve) to 
our simulations for the case of z = 1 
and z = 2. We find that although AGN 
sources preferentially lie on the AGN- 
sequence defined by Ivison et al. based on 
Mrk231, it is more acurately described as 
an extinction sequence at high (-2^2) red- 
shifts. PAH-dominated and continuum- 
dominated (i.e. AGN) sources are found 
on this sequence with roughly equal prob- 
ability, provided relatively high extinction 
levels (e.g. Ty ~ 10) are applied to both. 
We also note that the use of this plot for 
even rough redshift estimates is compli- 
cated by the large variations arising as the 
PAH 7.7 //m feature complex enters and 
leaves the 24 fira filter, meaning that high 
values of S'24/ 5*8 do not necessarily indicate 
low redshift. 

7.2. Future work 

In this work we have made use of the 
best mid-IR data available to date over 
a wide range of galaxy types. However, 
the Infrared Spectrograph (IRS) on Spitzer 
will soon improve on these, both in data 
quality and sample size. IRS does not ob- 
serve the A < 5 /im range (i.e. it will miss 
most of the direct stellar emission), but it 
does extend to 38 /im, which will allow us 
to include additional components to our 
model such as warm dust emission. We 
plan to make use of these future Spitzer 

^As we are only concerned with A > 1 /um, modest 
extinction cannot be excluded. 
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data in order to repeat this analysis in a 
complementary manner. 

Further, we plan to extend this ap- 
proach to modeling galaxy SEDs to longer 
wavelengths in particular far-lR/sub-mm. 
The far-IR/sub-mm is, by comparison with 
the mid-IR, relatively simple and can be 
modeled empirically with only 1-2 param- 
eters as shown for the ISO key project 
galaxies by Dale ct al. (2001). The ob- 
served correlations between our model pa- 
rameters and the IRAS colors suggest that 
such an extension is feasible^. 

The full potential of the MCMC ap- 
proach has not yet been exploited. In par- 
ticular, work is in progress to assess to 
what degree we can estimate the N(z), lu- 
minosity functions and so on of observed 
distributions, which is only hinted at in 
Fig. 11. For such a study, a more complete 
investigation of the effect of simulated an- 
nealing on the resulting distribution as well 
as potential biases and numerical optimiza- 
tions will need to be addressed. 

Lastly, we would like to understand fur- 
ther how this parameterization of galaxies' 
SEDs relates to physical properties such 
as star-formation rate, dust properties, 
and metallicity. We plan to approach this 
primarily through detailed UV/optical- 
to-MIR SED modeling of local, resolved 
galaxies. 

This work uses ISO archive data, specif- 
ically galaxy PHT-S spectra reduced and 
graciously provided to us by N. Lu and 
M. Hur of the Infrared Processing and 

^An additional indication that such an extension 
arises naturally is the observation that the PAH 
and sub- mm emission of resolved galaxies correlate 
fairly strongly (Haas, Klaas, & Bianchi 2002). 



Analysis Center (IPAC). We are also grate- 
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Fig. 1. — Here we illustrate the idea be- 
hind our model, where the 'activity level' 
of the galaxies increases toward the top. 
At the bottom is NGC4374 overlaid with 
the 10 Gyrs old stellar population template 
plus a cool continuum (a ~ 3). In the 
middle is M51 overlaid with the same as 
the above only with warmer continuum 
(a ~ 2), and our PAH emission template. 
At the top is IC4329a, which is essentially 
just a pure power law (a ~ 1). Note that 
substructure in the spectra remains; how- 
ever, these components govern the bulk of 
the MIR colour variations. 
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Fig. 2. — A sub-sample of our galaxies used 
to illustrate the continuum of PAH/stars 
relative strength. The solid line shows the 
best-fit model, while the dashed line shows 
the power-law component separately. The 
right-hand column shows the accepted tri- 
als for the stars (x-axis) and PAH (?/-axis) 
amplitude parameters (each is divided by 
TivFj^ for ease of comparison) for each 
galaxy. The scale is the same for all. 
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Table 1 
Basic data for our sample 



Name 


A f Tonnnl 

x\,A[JZUUUJ 


uec [Jzuuuj 


\jldiSS 


z 


iOg(,D60/'-'100j 


Arp220 


15:34:57.28 


23:30:11.3 


LINER, HIT, Sy2 


0.01813 


-0.03 


CenA 


13:25:27.50 


-43:01:10.9 


SO, Sy2 


0.00183 


-0.35 


ES0317-G23 


10:24:42.50 


-39:18:21.1 


SBa 


0.00965 


-0.28 


IC860 


13:15:03.53 


24:37:07.9 


Sa, HII 


0.01116 


-0.00 


IC883 


13:20:35.34 


34:08:22.0 


Im, Hll, LINER 


0.02335 


-0.21 


IC1459 


22:57:10.46 


-36:27:44.5 


E3, LINER 


0.00564 


-0.36 


IC4329a 


13:49:19.15 


-30:18:34.24 


SAO+, Syl.2 


0.01605 


0.09 


IC4595 


16:20:44.26 


-70:08:35.2 


Sc 


0.01137 


-0.42 


M51 


13:29:52.66 


47:11:42.6 


SAbc, HII, Sy2.5 


0.00154 


-0.63 


M81 


09:55:33.23 


69:03:55.0 


SAab, LINER, Syl.8 


0.000133 


-0.67 


M83 


13:37:00.46 


-29:51:55.7 


SABc, HII, Sbrst 


0.00172 


-0.33 


M87 


12:30:49.51 


12:23:28.6 


E, NLRG Sy 


0.00436 


-0.02 


Mrk231 


12:56:14:18 


56:52:25.26 


SAc, Syl 


0.04217 


0.02 


Mrk331 


23:51:26.24 


20:35:08.2 


S?, HII, Sy2 


0.01848 


-0.06 


Mrk359 


01:27:32.2 


19:10:39.5 


SBOa, Syl.5 


0.01739 


-0.19 


Mrk520 


22:00:41.43 


10:33:07.5 


Sb, Syl.9 


0.02661 


-0.21 


Mrk817 


14:36:22.08 


58:47:39.0 


SBc, Syl.5 


0.03145 


-0.03 


NGC0701 


01:51:03.71 


-09:42:10.2 


SBc, Sbrst 


0.00610 


-0.34 


NGC0863 


02:14:33.54 


00:46:00.4 


SAa, Syl.2 


0.02638 


-0.48 


NGC1022 


02:38:32.47 


-06:40:39.0 


SBa, HII, Sbrst 


0.00485 


-0.13 


NGC1056 


02:42:48.34 


28:34:26.8 


Sa, HII, Sy 


0.00515 


-0.25 


NGC1068 


02:42:40.83 


-00:00:47.0 


SAb, Syl, Sy2 


0.00379 


-0.10 


NGC1097 


02:46:19.06 


-30:16:28.0 


SBb, Syl 


0.00425 


-0.28 


NGC1125 


02:51:39.51 


-16:39:08.4 


SAB, Sy2 


0.03175 


-0.07 


NGC1222 


03:08:56.81 


-02:57:17.6 


SO, HII 


0.00818 


-0.06 


NGC1326 


03:23:56.40 


-36:27:49.6 


SBO/a, LINER 


0.00454 


-0.25 


NGC1377 


03:36:39.01 


-20:54:06.2 


SO, HH 


0.00598 


0.10 


NGC1385 


03:37:28.22 


-24:30:04.0 


SBcd 


0.00498 


-0.33 


NGC1546 


04:14:36.60 


-56:03:38.9 


SAO+? 


0.00426 


-0.58 


NGC1569 


04:28:11.14 


64:50:52.25 


IBm, Sbrst, Syl 


0.000244 


-0.02 


NGC3227 


10:23:30.60 


19:51:55.2 


SAB, Syl.5 


0.00386 


-0.35 


NGC3705 


11:30:05.82 


09:16:36.3 


SABab, LINER, HII 


0.00340 


-0.46 


NGC3783 


11:39:01.78 


-37:44:19.6 


SBa, Syl 


0.00973 


-0.18 


NGC3949 


11:53:41.42 


47:51:31.7 


SAbc, HII 


0.00267 


-0.38 


NGC4027 


11:59:30.59 


-19:15:48.1 


SBdm 


0.005574 


-0.39 


NGC4051 


12:03:09.58 


44:31:52.9 


SABbc, Syl.5 


0.00234 


-0.53 


NGC4102 


12:06:23.37 


52:42:41.2 


SABb?, HII, LINER 


0.00282 


-0.17 


NGC4194 


12:14:10.10 


54:31:41.9 


IBm, BCG, HII 


0.00836 


-0.08 


NGC4253 


12:18:26.43 


29:48:46.0 


SBa, Syl.5 


0.01293 


-0.06 


NGC4365 


12:24:28.87 


07:19:04.6 


E3 


0.00415 


-1.17 
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Fig. 3. — The histograms of the IRAC 
fluxes obtained from the residuals. Here 

the dashed hne indicates the scatter ex- 
pected from the uncertainties in the data, 
while the dotted line is the Gaussian cor- 
responding to the mean and a of the his- 
tograms. The positive excess noticeable in 
the 4.5 fim plot may be due to the contri- 
bution of various strong lines which emit 
there (e.g. Bra, CO), although some resid- 
ual poor calibration can also be the cause, 
as suggested by the presence of a nega- 
tive excess as well. The quality of the 
data does not allow us to distinguish these. 
The excess scatter in the four bands is: 
0-3.6 = 10 mJy, £74.5 = 8 mJy, (75.8 = 18 mJy, 
and (jg.o = 9 mJy. 
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Fig. 4. — The marginalized distributions 
of the continuum slope a. The top six 
panels show these for a representative sub- 
sample, while the bottom panel shows the 
total distribution for the three spectral 
types, where the solid line is the stellar- 
dominated distribution, the dotted line is 
for PAH-dominated, and the dot-dashed 
line is for continuum-dominated. Note the 
clear difference in the peaks as well as the 
typically double-peaked distributions. See 
Section 5.2 for discussion. 
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Fig. 5. — Here PAH/continuum/stars 
mean the best-fit components integrated 
over the 3-12 //m range. Only z > 0.004 
galaxies are used, to minimize aperture 
bias. Stars mark three SF-dominated 
active galaxies, while circles show three 
AGN-dominated active galaxies. See Sec- 
tion 5.2 for details. 



Fig. 6. — Here we show the Monte 
Carlo generated IRAC color-color diagram, 
where red indicates starlight-dominated, 
green is PAH-dominated, and blue desig- 
nates continuum-dominated. The tracks 
show the effect of AGB dust shells on the 
IRAC colors for a range of stellar popula- 
tion ages (given in log[Myr]). 
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Fig. 7.— The SSP spectra (Piovan, Tan- 
talo, Sz Chiosi 2003) corresponding to the 
age track in Fig. 6 overlaid with represen- 
tative galaxy spectra. See Section 6.1 for 
discussion. 
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Fig. 8.— The redshift evolu- 

tion of three model spectra representa- 
tive of the three major types of mid- 
IR spectra: continuum-dominated, PAH- 
dominated; and stellar-dominated. We fo- 
cus on the log(5'8.o/'S'4.5) vs. log(5'5.8/-S'3.6) 
projection in this work, since it appar- 
ently provides the best separation be- 
tween the PAH-dominated and continuum- 
dominated sources in particular. 



Fig. 10. — The z = 0-2 color evolution 
based on our simulations. The color 
scheme is as in Fig. 6. The top panels 
are divided into regions which preferen- 
tially select continuum-dominated sources 
at z — 2 (region 1), PAH-dominated 
sources at 2; ~ 0.05 — 0.3 (region 2), stellar- 
and PAH-dominated sources at 2;~0.3- 
1.6 (region 3), and same at 2; ^1.6 (region 
4). The bottom panels show simulations in 
the range 2; = 1-2 only. The z = region 
marks the spread of colors for our sample, 
while the 2; < 1 is our best guess based on 
the assumption that the lg(5'24/S'5,8) will 
not drop below this range before the 24 fim 
filter begins to sample the 10 /^m dip at 
zr^l.5. Stellar-dominated systems here 
are off-scale (see Fig. 9). See Section 6.4 
for a discussion. 



Fig. 9. — The color evolution of the 
stellar-dominated, continuum-dominated, 
and PAH-dominated sources as a function 
of redshift. Note in particular the influence 
of the 7.7 iim PAH feature a,t z ^ 0.5, and 
the 1.6 /xm peak at z > 1.3. The scatter 
represents ± 2 cr. See Section 6.2 for fur- 
ther details. 
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Fig. 11.— top: The FLS IRAC color- 
color plot, bottom: A simulation from our 
model. See Section 7.1 for details and dis- 
cussion. 



Fig. 12. — The effect of extinction. Here 
as usual green indicates PAH-dominated 
sources and blue is continuum-dominated 
sources. The red filled circles are the 
MAMBO sources from Ivison et al. (2004). 
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Table 1 Cnnfmiicd 



Name 


RA[J2000] 


Dec[J2000] 


Class"* 


z 


log(S'6o/S'ioo) 


NGC4374 


12:25:03.09 


12:53:10.5 


El, LERG, LINER 


0.00354 


0.37 


NGC4418 


12:26:54.65 


00:52:40.5 


SABa, Sy2 


0.00727 


0.09 


NGC4490 


12:30:36.88 


41:38:23.4 


SBd 


0.00188 


0.25 


NGC4507 


12:35:36.67 


39:54:33.16 


SABab, Sy2 


0.01180 


0.10 


NGC4569 


12:36:49.78 


13:09:46.2 


SABab, LINER, Sy 


0.002148 


0.50 


NGC4593 


12:39:39.35 


05:20:38.8 


SBb, Syl 


0.00900 


- 0.29 


NGC4691 


12:48:13.47 


03:19:59.1 


SBO/a, HII 


0.00137 


-0.15 


NGC5253 


13:39:55.75 


-31:38:30.8 


Im, HII, Sbrst 


0.00135 


-0.011 


NGC5433 


14:02:36.03 


32:30:37.5 


Sdm 


0.01452 


-0.25 


NGC5548 


14:17:59.48 


25:08:12.11 


SAO/a, Syl. 5 


0.01717 


-0.18 


NGC5866 


15:06:29.41 


55:45:47.2 


SO, HII/LINER 


0.00224 


-0.57 


NGC5915 


15:21:33.16 


-13:05:32.5 


SBab 


0.00758 


-0.17 


NGC5962 


15:36:31.70 


16:36:31.8 


SAc, HH 


0.00653 


-0.42 


NGC6000 


15:49:49.28 


29:23:13.1 


SBbc, HII 


0.00732 


-0.17 


NGC6753 


19:11:23.29 


-57:02:55.4 


SAb 


0.01042 


-0.46 


NGC7218 


22:10:11.71 


-16:39:35.7 


SBc 


0.00554 


-0.34 


NGC7418 


22:56:35.90 


-37:01:45.7 


SABcd 


0.00482 


-0.54 


NGC7469 


23:03:15.58 


08:52:25.9 


SABa, Syl. 2 


0.01632 


-0.13 


NGC7592 


23:18:22.10 


-04:24:59.8 


LINER, Sy2 


0.02444 


-0.12 
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